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INTRODUCTION 


The purpose of this paper is to present an easily auto- 
mated technique for computing a continuous function to approx 
imate a discrete set of data. An algorithm is out?.ined 
describing a weighted least squares solution involving the 
approximation for the dependent variable an- the desired 
derivative of the approximating function. This pape*' will 
be concerned indirectly with numerical differentiation. The 
algorithm can be applied directly but requires a few more 
computations than standard least squares. 

NUMERICAL APPROXIMATIONS 


Discrete data sets are often replaced by functional 

approximations for the purpose of interpolation, integration, 

* 

or differentiation. None of these processes are precise. 

For example, given any finite set of ordered pairs of real 
numbers » i = 1,N , and any finite interval 

[a,b] , there is an infinite collection of continuous, 
differentiable, integrable functions F defined on [a,b] 
such that, for every function f in F , f(x^) = y i . 

Moreover, for every x Q not in the original data set x^^ , 
there is an infinite collection of functions in F all 


having distinct values at x Q , distinct derivatives at 
x Q , and distinct integrals over [a,b] . This implies that 
any numerical algorithm for interpolation, differentiation, 
or integration of finite data sets yields an answer that is 
at best a "correct answer" among an infinite number of 
"correct answers". Because of this fact, most texts describe 
algorithms that require some additional properties of the 
data or present algorithms which appeal to the intuition of 
the user. This paper utilizes the latter alternative. 


A MINIMIZATION TECHNIQUE 


In all that follows it will be assumed that a fixed set 

of data (x^y^) for i = 1,N has been defined. For a 

fixed set of functions >**2 » * * * *k an< ^ the ^ unct ^ on 
k 


r -E a i f i 


the least squares norm of f - y will be 


1=1 


defined by 


Ilf - y| I 



(i) 


and the derivative of f with respect to x will be denoted 
by f* . Least squares approximations of this type often 


exhibit "undesirable behavior" between the input data points. 


Lagrange and Hermite interpolating polynomials provide approx- 
imations exact at the input points, but they do not provide 
any smoothing of the data. The following algorithm is similar 
to the Hermite interpolation and the least squares smoothing 
process . 

Given two nonnegative weight functions wj(x) and 
w ? (x) , minimize 

I | w x (f - y ) + w 2 (f ' - y ) | | (2) 

A 

where y^ is the desired derivative of f at x^ . A 
function f minimizing (2) thus combines the requirements 
that f approximates y and the slope of the curve f(x) 
approximates the desired slope at each of the initial datum 
points . 

Rewriting (2) in the form 

| |wi f + w 2 f* - (wj y + w 2 y ) | | (3) 

and considering first the restriction of f to a polynomial 
of degree k-1 , we see that the minimization problem can 
easily be solved as follows: 

Define an Nxk matrix B where the last k-1 columns are 


defined by 


b i 1 B w i^ x i^ x ‘]” 1 + U-l) w 2 (x a ) xj" 2 (*0 

and the first column is defined 

b n e w i^ x i) (5) 

Let z be a column vector such that 

z i c w i^ x i^ y i + W 2^ X i^ y i 

T 

If B B is nonsingular let 

a «= (B T B)- 1 B T z (7) 

and if B* 1 B is singular let 

a « B + z (8) 

where B + is the generalized inverse of B . Then a func- 
tion f minimizing (3) is given by 

i ■, -- ■ 

i=l 


f ( x ) 


The only restriction needed for the class of functions 
to insure a simple solution is that for each function f A , 

f^ - fj for some ,1 . For examnle if f Is restricted to 

f(x) * a, + a„ x + a. x 2 + a, sin(2x) + a e cos(2x) 

then the first three columns of B are defined as in 00 

and (S), ond the other columns are defined by 

b ji| B w i( x i) sin(2x 1 ) + ?w i> (x 1 ) cos(2x 1 ) 
b il> B w i^ x i^ cos(2x i ) - 2w 2 ( x i ) sin(2x A ) 

The vector z remains the same and the solution is obtained 
in exactly the same manner. 

CONCLUS ION 

This algorithm has been programed and used for several 
applications. It should be pointed out that orthogonal poly- 
nomials can be used to obtain the solution to Equation (7) 
v/hen a polynomial approximation is computed for Wj ■ w 2 ■ 1 . 
Also, notice that when Wj (or w 2 ) is constant valued only 


c 


one value should be reserved in computer memory. Another 
option that can be implemented is to use the average value 
of the slopes of the left and right straight line segments 
between adjacent points for the value of . In this case, 

core storage for the y array need not be reserved. 


